Quick exercise in plotting UK national crime database data using OpenStreetMap, and fitting a contour (kernel density estimator) function for 2-D spatial smoothing. Will end up looking something like:
Example of crime data plotted onto OpenStreetMap, with plotting function
# First the required packages.
# The top-level ones are OpenStreetMap, MASS and ks
# Other dependencies below this include:
# X11, rgl, rgdal
# On linux this will mean something like:
# sudo apt-get update && sudo apt-get install libcgal-dev libglu1-mesa-dev libglu1-mesa-dev libx11-dev gdal-bin gdal-config
# I googled the shit out of this, BTW:
# https://stackoverflow.com/questions/15248815/rgdal-package-installation
# install.packages("OpenStreetMap")
# install.packages("ks")
# install.packages("MASS")
# load required libs
library(OpenStreetMap)
library(MASS)
library(ks)
# Get the crime data itself
# From e.g. https://wwww.police.uk/surrey/RUCM/crime/stats
# Saved to ~/Downloads/c1b1d294-4b4b-4a47-8f9c-4dcd7af9e56a.csv
crime <- read.csv("~/Downloads/c1b1d294-4b4b-4a47-8f9c-4dcd7af9e56a.csv")
head(crime)
## Crime.ID Month
## 1 554fd7ba6fc67d388c02a1be61837a2c819843957f1ed0852a705647e0b2e14e 2017-11
## 2 554fd7ba6fc67d388c02a1be61837a2c819843957f1ed0852a705647e0b2e14e 2017-11
## 3 554fd7ba6fc67d388c02a1be61837a2c819843957f1ed0852a705647e0b2e14e 2017-11
## 4 554fd7ba6fc67d388c02a1be61837a2c819843957f1ed0852a705647e0b2e14e 2017-11
## 5 554fd7ba6fc67d388c02a1be61837a2c819843957f1ed0852a705647e0b2e14e 2017-11
## 6 554fd7ba6fc67d388c02a1be61837a2c819843957f1ed0852a705647e0b2e14e 2017-11
## Latitude Longitude Location Category
## 1 51.38579 -0.489942 On or near Meadow View Anti-social behaviour
## 2 51.38944 -0.511582 On or near Pyrcroft Road Anti-social behaviour
## 3 51.39203 -0.503233 On or near London Street Anti-social behaviour
## 4 51.38970 -0.515498 On or near Vincent Road Anti-social behaviour
## 5 51.39067 -0.504155 On or near Beomonds Row Anti-social behaviour
## 6 51.39601 -0.510032 On or near Cowper Close Anti-social behaviour
## Outcome.status
## 1
## 2
## 3
## 4
## 5
## 6
The timestamps here aren’t important but might as well convert them:
crime$datmon <- as.Date(paste(crime$Month,'-01',sep=''))
Now we’re ready to plot the map itself. Notes from https://www.r-bloggers.com/plot-maps-like-a-boss/.
# get the map
chertsey <- (openmap(c(51.40,-0.52), c(51.37,-0.48)))
plot(chertsey)
Looks cool!
Note, map() takes a top-left then bottom-right pair, e.g. map(c(x1,y1),c(x2,y2)) where x will be longitude and y latitude co-ords. This produces a plotting device you can add to like any other plot() call, BUT without transformation/projection the points returned in e.g. OpenStreetMap will be in (I think) metres N and E from the equator/meridien (I think). This means that to plot points on them directly you have to do some guessing:
plot(chertsey)
abline(v=seq(-57886,0,length.out=30),col="green",lwd=4)
# wtf?
The lines appear at seemingly crazy intervals (“aren’t these lat/long co-ordinates? What does a longitude of ‘-57,000’ mean?!”). This is because, importantly, the crime dataset lat/long points refer to a flat 2D surface, whereas we need to project to the OpenStreetMap spheroid.
Finally we can plot and visualise our surface:
#transform map from OpenStreetMap co-ords to 2D lat-long projection
map_ll <- openproj(chertsey,projection="+proj=longlat")
# plot map
plot(map_ll)
# plot individual crimes on the map in big fat purple
points(Latitude~Longitude,data=crime,col="purple",pch=19,cex=1.5)
# create a contour (surface) plot using KDE
crm <- kde2d(crime$Longitude,crime$Latitude)
# add the contour
contour(crm,add=TRUE,col=rev(heat.colors(11)),lwd=5)
Hey presto!
BTW we might also wonder if particular categories of crime occur in particular places. I can’t be arsed to do this in full but clearly we’d subset the data, etc. Quick check with boxplot() tho:
boxplot(Latitude~Category,data=crime)
boxplot(Longitude~Category,data=crime)
So… there aren’t really any particular locations more/less prone to any type of crime.